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Предложен подход к моделированию изменений земной поверхности, основанный на вероятностном формировании правил 
функционирования клеточных автоматов. Приводятся результаты исследований предлагаемого подхода, проведенные на мо- 
дельных разновременных изображениях. Обсуждаются результаты апробации этого подхода в случае анализа динамики зем- 
ной поверхности фрагмента территории Ханты-Мансийского автономного округа поданным космосъемки. 


Введение 

Земная поверхность представляет собой слож- 
ную систему, а моделирование ее изменений явля- 
ется сложным процессом, на который влияет боль- 
шое количество всевозможных факторов. Для того, 
чтобы представить некоторый интересующий 
фрагмент земной поверхности в качестве объекта 
моделирования, необходимо получить изображе- 
ние этого фрагмента в виде матрицы, каждый 
элемент которой соответствует конкретному типу 
земной поверхности. Фактически эта матрица 
представляет собой тематическую карту такого 
фрагмента земной поверхности. Такие тематиче- 
ские карты, необходимые для моделирования из- 
менений земной поверхности, часто получают с 
помощью методов дистанционного зондирования 
и систем автоматизированной интерпретации аэ- 
рокосмических изображений (АИ). 

На сегодняшний день одним из наиболее эф- 
фективных и широко используемых способов мо- 
делирования изменений земной поверхности, учи- 
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тывающих пространственное взаимодействие меж- 
ду элементами подобного изображения, является 
применение клеточных автоматов (КА) [1]. Одним 
из ключевых факторов, влияющих на эффектив- 
ность их использования, является определение в 
каждом конкретном случае правил взаимодействия 
элементов КА или правил функционирования КА. 

В данной работе предлагается подход, при кото- 
ром правила функционирования КА имеют более 
сложную вероятностную природу, чем при тради- 
ционном детерминистском подходе. Предлага- 
емый подход учитывает особенности моделируемо- 
го изображения земной поверхности и ведет к по- 
лучению более эффективных правил функциони- 
рования для каждого используемого при моделиро- 
вании КА. Все это позволяет проводить моделиро- 
вание изменений земной поверхности более эф- 
фективно и получать более точные результаты мо- 
делирования по сравнению с результатами, полу- 
ченными в рамках существующих моделей земной 
поверхности. 
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Теоретические основы моделирования изменений 

земной поверхности 

Один из наиболее перспективных и рациональ- 
ных способов моделирования изменений земной 
поверхности основан на использовании статистиче- 
ского аппарата марковских цепей [1]. В этом случае 
вероятностная информация об изменениях земной 
поверхности содержится в матрице вероятностей 
переходов Р=[д,], каждый элемент которой является 
вероятностью перехода типа земной поверхности со, 
в тип со } , где ц = 1 М— число типов земной по- 
верхности на изучаемой территории, выявленное на 
этапе интерпретации АИ. Количественная инфор- 
мация об изменениях содержится в матрице факти- 
ческих переходов территорий М факг =[т|’] и матрице 
ожидаемых переходов территорий М ожид =[т“]. Ма- 
трица М факт содержит количество элементов изобра- 
жения (пикселей), перешедших из типа земной по- 
верхности а>і в тип со р а М 0Ж1ІД - количество пикселей 
изображения, которые должны поменять свое со- 
стояние (значение) в результате процесса моделиро- 
вания изменений земной поверхности. 

Множество КА, с использованием которых осу- 
ществляется моделирование, формируется при после- 
довательном сканировании исходного изображения. 


При этом каждый КА можно представить в виде ква- 
дратной матрицы М СА =[с^] порядка д. Значение цен- 
трального элемента с' № каждой такой матрицы зависит 
некоторым образом от значений всех элементов этой 
матрицы с' и =/ СА (с 1Ь с 12 ,...,с и ,...,с й ). Полученное при 
этом значение с' И , является элементом результирующе- 
го изображения. При традиционном подходе поиск 
функции / СА чаще всего основан на использовании ка- 
кого-либо простого правила (например, используется 
принцип мажоритарности), найденная функция ис- 
пользуется для всех КА моделируемого изображения. 
В то же время в предлагаемом подходе к моделирова- 
нию изменений земной поверхности применяется та- 
кой способ формирования функции^, при котором 
для каждого КА используется три различных вероят- 
ностных компонента ( р р^, д/ ірст ), что позволяет 
учесть особенности формирования правил взаимодей- 
ствия каждого используемого КА. При этом каждый 
из компонентов представляет собой соответствующие 
вероятности р ц перехода типа земной поверхности ш, в 
тип Юр Рассмотрим более подробно особенности рас- 
чета каждой такой вероятности. 

Компонент, учитывающий факторы влияния. В 
настоящее время при моделировании изменений 
заселенных территорий (с помощью моделей, из- 
вестных как ІІгЬап ОгоѵіН Модек) широко исполь- 
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зуется подход, основанный на применении так на- 
зываемых карт соответствия [1]. Эти карты несут 
информацию о вероятности заселения каждой точ- 
ки земной поверхности изучаемого фрагмента тер- 
ритории в зависимости от таких факторов влияния 
как расстояние до автомагистралей, расстояние до 
крупных торговых центров, расстояние до источ- 
ников питьевой воды, величина уклона местности 
и т.п. Для применения подобного подхода при мо- 
делировании изменений земной поверхности 
необходимо построение таких карт соответствия 
для всех типов земной поверхности, представлен- 
ных на тематической карте. Учитывая возможные 
сложные взаимозависимости между подобными 
факторами, задача построения подобных карт со- 
ответствия для каждого типа земной поверхности 
представляется сложной и трудно формализуемой. 

В настоящее время одним из наиболее эффективных 
способов, позволяющих использовать подобные нефор- 
мализованные исходные данные при моделировании 
изменений, является применение искусственных ней- 
ронных сетей (ИНС). Наиболее распространенной то- 
пологией для решения подобной задачи является много- 
слойный персептрон (МСП) с одним скрытым слоем. 
Применим технологию использования ИНС, аналогич- 
ную той, что описана в [2], которая основана на нейро- 
сетевых КА. Суть этой технологии заключается в том, 
что для моделирования роста заселенных территорий 
применяются КА, правила функционирования которых 
определяются с помощью ИНС. Особенности примене- 
ния МСП в рамках предлагаемого подхода при обуче- 
нии и получении отклика изображены на рис. 1. 

В общем случае применяемый МСП имеет 
К=М Р входных и М выходных нейронов, где М - 
число типов земной поверхности на разновремен- 
ных тематических картах, Г - число факторов 
влияния. Процесс получения необходимой вероят- 
ностной информации с помощью нейросети, учи- 
тывающей выявленные факторы влияния, заклю- 
чается в том, что каждый аксон /-ого выходного 
нейрона интерпретируется как вероятность пере- 
хода анализируемого пикселя типа со , исходного 
изображения в тип щ (рІІ' жт ). 

Компонент, упитывающий вероятностную соста- 
вляющую. Как упоминалось выше, привлечение веро- 
ятностной информации о процессах, происходящих 
на земной поверхности при моделировании часто ба- 
зируется на аппарате марковских случайных цепей. 
Одной из основополагающих составляющих этого 
аппарата является матрица вероятностей переходов 
Р =[д,], традиционно получаемая на основе двух раз- 
новременных тематических карт изучаемого фраг- 
мента земной поверхности. Таким образом, вычисле- 
ние второго компонента в рамках предлагаемого под- 
хода основывается на элементах матрицы вероятно- 
стей переходов так, что в некоторой окрестности ве- 
роятность перехода типа со, земной поверхности в тип 
о, зависит не только от вероятности р,,, но и от коли- 
чества элементов типа щ в этой окрестности. Для 
каждого типа земной поверхности в анализируемой 


окрестности, совпадающей с элементами текущего 
КА, вычисляется вероятность рР'=п-р^ у-1,2,...,/я, 
где т - количество типов земной поверхности в ана- 
лизируемой окрестности, а количество элементов 
типа оу в анализируемой окрестности. 

Компонент, учитывающий пространственную со- 
ставляющую. Третий компонент формирования пра- 
вил функционирования КА базируется на использо- 
вании пространственных характеристик каждого из 
типов земной поверхности, основанных на вычисле- 
нии фактора насыщенности, предложенного в [3]. 
Использование факторов насыщенности происходит 
следующим образом. Для каждого со, формируется 
соответствующий вектор Г," ас ={/' і=1,...,М, где 
М- количество типов земной поверхности. Этот век- 
тор Р™ несет информацию о насыщенности того или 
иного типа земной поверхности на всем изображе- 
нии в целом. Затем, для каждой точки изображения 
рассчитывается локальный вектор насыщенности 
Р™, после чего вероятность р і , прст определяется как 
р^~ 1/г/(Р ; нас ,Б’, жж ), где й - евклидово расстояние меж- 
ду двумя этими векторами насыщенности. 

После определения каждого из вероятностных 
компонентов, результирующую вероятность пере- 
хода каждого из элементов исходного изображения 
в выбранной окрестности текущего КА можно 
представить как р^^^^р^+р^+р,^ ■ 

Немаловажной необходимой деталью при про- 
ведении моделирования изменений земной по- 
верхности является не только определение вероят- 
ности перехода р*“", но и порядок изменения эл- 
ементов результирующего изображения. Этот по- 
рядок в рамках предлагаемого подхода определяет- 
ся следующим образом. Самыми первыми изменя- 
ются элементы изображения, результирующая ве- 
роятность для которых р^ ш имеет наибольшее зна- 
чение, а самыми последними - те элементы изо- 
бражения, для которых эта вероятность наимень- 
шая. Предварительные результаты исследования 
эффективности применения такого ранжирования 
показывают, что оно позволяет получать значи- 
тельно более адекватные результаты моделирова- 
ния изменений земной поверхности [4]. 

Исследование предлагаемого подхода с использова- 
нием модельных разновременных изображений 

Исследование эффективности предлагаемого 
подхода с точки зрения точности было проведено с 
использованием модельных разновременных тема- 
тических карт. Модельные данные представляли 
собой 30 различных комбинаций из трех модель- 
ных изображений, имитирующих разновременные 
АИ, на которых изменение границ и площадей ти- 
пов земной поверхности происходит аналогично 
изменениям, наблюдаемым на реальной земной 
поверхности. Моделирование осуществляется по 
двум изображениям, а точность оценивается по 
третьему изображению с помощью матрицы пере- 
путывания и интегрального критерия точности - 
каппа индекса согласия (КИС) [4]. Для получения 
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статистически достоверных результатов экспери- 
менты тридцатикратно повторялись, а в качестве 
оценки КИС результата эксперимента принима- 
лось его среднее значение. 



1 , 2->3 2 , 3->4 4 , 3->2 3 , 2->1 

Пары разновременных изображений 


Рис. 2. Точность прогнозных карт, полученных с использова - 
нием различных составляющих в правилах функцио- 
нирования КА: 1) Ру 2) р ® р , 3 ) р“ р и р , " рсг , 4) точность 
карты, полученной с помощью модуля СА_Магкоѵ 


В отличие от вероятностных характеристик р ^ 
и рГ, рассчитываемых исключительно по разно- 
временным тематическим картам, вероятность рА ш 
определяется с использованием дополнительной 
априорной информации о процессах, происходя- 
щих на исследуемом фрагменте земной поверхно- 
сти. Для простоты, исследования проводились при 
условии отсутствия влияния априорной информа- 
ции, поэтому вероятность рт ш не учитывалась. От- 
метим также, что проводились сравнительные ис- 
следования эффективности построения прогноз- 
ных карт, выполненные с помощью существующего 
модуля моделирования изменений земной поверх- 
ности СА_Магкоѵ, представленного в растровой 
геоинформационной системе (ГИС) МгІ8І32 (США) 
[7]. На рис. 2 в качестве примера представлена часть 
результатов проведенных исследований. На оси 
абсцисс показаны различные виды комбинаций че- 
тырех модельных изображений, имитирующих раз- 
новременные АН. Так, запись 1 ,2 — >3 означает, что 
расчет параметров и построение прогнозной карты 
проводилось на основе модельных изображений 
№ 1 и 2, а точность полученной прогнозной карты 
оценивалась по изображению № 3. 

Результаты проведенных исследований показа- 
ли высокую значимость вероятностных характери- 
стик р’ с " и р"^ при формировании правил функци- 
онирования КА, а также то, что созданное про- 
граммное обеспечение (ПО), реализующее предло- 
женный подход, позволяет строить значительно 
более точные прогнозные карты, чем модуль 
СА Магкоѵ ГИС Іс1гі$і32. 


Применение предлагаемого подхода для анализа 
динамики земной поверхности фрагмента террито- 
рии Ханты-Мансийского автономного округа 

В ГИС лаборатории ТПУ были программно ре- 
ализованы предлагаемый подход к моделированию 
изменений земной поверхности и соответствую- 
щие алгоритмы, а также необходимые функции 
предварительной обработки и автоматизированной 
интерпретации АИ [8,9]. Апробация предлагаемого 
подхода проводилась при решении практически 


важной задачи анализа динамики земной поверх- 
ности фрагмента территории Ханты-Мансийского 
автономного округа (ХМАО). Решение задачи по- 
зволит выявить тенденции в изменении земной по- 
верхности в зоне магистрального нефтепровода, 
проходящего через исследуемую территорию. 

Для решения указанной задачи были использова- 
ны снимки системы ДЭЗ Ьапсікаі (ЕТМ+), получен- 
ные Центром приема и обработки данных ДЗЗ Югор- 
ского НИИ информационных технологий на 10 авгу- 
ста 1999 г., 6 сентября 2000 г. и 3 июля 2002 г. Все сним- 
ки имели уровень обработки 1С и были переведены в 
проекцию ЦТМ на эллипсоиде \У0884 по данным ор- 
битальной привязки. Учитывая точность проведен- 
ной геометрической коррекции и разрешение исход- 
ных АИ, рабочий масштаб тематических карт будет 
соответствовать картам масштаба 1 :200 000. 

Для анализа эффективности предлагаемого 
подхода и созданного ПО в целом в качестве базо- 
вых АИ для проведения интерпретации и модели- 
рования использовались только АИ 1999 и 2000 гг. 
В этом случае АИ 2002 г. можно принимать в каче- 
стве тестового (эталонного) для оценки качества 
результатов интерпретации и моделирования в ви- 
де прогнозной карты 2002 г. 

В итоге для решения задачи анализа динамики 
земной поверхности фрагмента территории ХМАО 
по имеющимся разновременным АИ была предло- 
жена технология получения прогнозных карт, ос- 
нованная на использовании разработанного ПО. 
Рассмотрим основные этапы этой технологии, ко- 
торая позволит дать представление о тенденциях, 
имеющих место на исследуемом фрагменте земной 
поверхности. 

Этап 1. Из имеющихся космических снимков, 
покрывающих территорию ХМАО, была составле- 
на композиция, которая позволила выявить такие 
фрагменты для всех трех космических снимков, ко- 
торые не закрыты облачным покровом и содержат 
общий фрагмент земной поверхности на исследуе- 
мой территории. Учитывая небольшие вариации 
космических снимков по фенологической фазе, 
имеется возможность построения обучающих вы- 
борок для одних и тех же типов земной поверхно- 
сти, представленных на изучаемой территории в 
основном различными классами растительности. 

Этап 2. Выбранные фрагменты космических 
снимков, которые являются общими для всех раз- 
новременных АИ, была подвергнуты неконтроли- 
руемой классификации по методу І80БАТА. Для 
этого использовалось ПО системы ЕК Маррег 5.5, а 
количество кластеров при неконтролируемой клас- 
сификации было экспериментально принято рав- 
ным десяти. С учетом результатов этой классифи- 
кации, а также имеющихся данных наземных на- 
блюдений, для каждого из космических снимков 
были построены обучающие выборки для таких ти- 
пов земной поверхности как пойменные участки, 
кустарники, мшистые болота и вода. Затем все кос- 
мические снимки были подвергнуты тематической 
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Рис. 3. Фрагменты разновременных тематических карт, полученных с помощью разработанного ПО: а) 1999 г., 6) 2000 г., в) 2002 г. 


интерпретации с помощью созданного ПО. Это по- 
зволило получить разновременные растровые те- 
матические карты, на которых представлены тема- 
тические классы, соответствующие построенным 
обучающим выборкам (рис. 3). 

Этап 3. Заключается в построении вероятност- 
ных карт соответствия с помощью модулей ПО ГИС 
Мгш32, встроенных в систему. Для этого необходима 
дополнительная априорная информация о вероят- 
ности замещения одного типа земной поверхности 
другим в каждой точке изучаемой территории {р~ т ), 
на основании которой с использованием функций 
пространственного анализа могут быть построены 
карты соответствия. В данном случае такая априор- 
ная информация отсутствовала и не учитывалась на 
дальнейших этапах предложенной технологии. 

Этап 4. После того как тематические карты по- 
строены, используя ПО подсистемы моделирова- 
ния изменений земной поверхности, рассчитыва- 
ются численные характеристики изменений, про- 
исходящих на земной поверхности. Для этого фор- 
мируются матрица вероятностей переходов (табл. 1) 
и матрица фактических переходов (табл. 2). 


Таблица 1. Матрица вероятностей переходов 



Пойменные 

участки 

Кустарники 

Вода 

Мшистые 

болота 

Пойменные 

участки 

0,375998 

0,221505 

0,001555 

0,400942 

Кустарники 

0,101386 

0,843558 

0,000232 

0,054824 

Вода 

0,094620 

0,005566 

0,812616 

0,087199 

Мшистые 

болота 

0,043285 

0,017094 

0,002255 

0,937366 


Расчет основывается на тематических картах, 
полученных с космических снимков 1999 и 2000 гг. 
Рассчитывается также матрица ожидаемых перехо- 
дов (табл. 3), необходимая для получения темати- 
ческой карты на следующий момент времени - 
200 1 г. Отметим, что матрица ожидаемых переходов 
и матрица фактических переходов формируются из 
пикселей соответствующих тематических карт 
(изображений), поэтому в табл. 2 и 3 приводятся 
элементы матрицы в пикселях. 


Таблица 2. Матрица фактических переходов 


Тип 

Пойменные 

участки 

Кустарники 

Вода 

Мшистые 

болота 

Пойменные 

участки 

130818 

77247 

554 

139822 

Кустарники 

44624 

375793 

128 

23614 

Вода 

1073 

58 

8560 

867 

Мшистые 

болота 

38144 

15356 

2046 

834497 


Используя информацию о вероятностных и 
других количественных тенденциях изменения ис- 
следуемой территории, содержащуюся в табл. 1-3, 
рассчитывается матрица средних значений насы- 
щенности, которая необходима для определения 
правил функционирования КА при моделирова- 
нии. Результаты этого расчета приведены в табл. 4. 

Этап 5. После того, как вся необходимая веро- 
ятностная и другая количественная информация 
получена, проводится непосредственно моделиро- 
вание изменений земной поверхности и получение 
прогнозной карты 2002 г. (рис. 4, б). В качестве про- 
межуточного результата моделирования, позво- 
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Рис. 4. Фрагменты прогнозных карт, полученные с помощью ПО системы: а) 2001 г, б) 2002 г. 
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Технические науки 


ляющего оценить интересующие нас изменения 
земной поверхности 2001 г., получена прогнозная 
тематическая карта 2001 г. (рис. 4, а). 


Таблица 3. Матрица ожидаемых переходов 


Тип 

Пойменные 

участки 

Кустарники 

Вода 

Мшистые 

болота 

Пойменные 

участки 

62924 

96774 

1000 

187743 

Кустарники 

55534 

328251 

337 

60038 

Вода 

1315 

348 

6944 

1949 

Мшистые 

болота 

51835 

35857 

3642 

798709 


Этап 6. Для того чтобы осуществить дополни- 
тельную (сравнительную) оценку полученной с по- 
мощью ПО системы прогнозной карты на 2002 г., бы- 
ла сформирована еще одна прогнозная карта на 
2002 г. Эта карта была получена с использованием 
аналогичных входных данных (разновременных те- 
матических карт), что и в первом случае, но с помо- 
щью модуля САМагкоѵ, представленного в ГИС Мгі- 
8І32. Также во втором случае в качестве дополнитель- 
ного, промежуточного результата моделирования бы- 
ла получена прогнозная карта на 2001 г. Заметим, что 
при получении прогнозных карт использовались па- 
раметры модуля СА_Магкоѵ, задаваемые по умолча- 
нию, а размер анализируемой окрестности КА был 
аналогичен размеру в первом случае и составлял 5x5. 


Таблица 4. Матрица значений насыщенности исходной тема- 
тической карты на 2000 г. 



Пойменные 

участки 

Кустарники 

Вода 

Мшистые 

болота 

Пойменные 

участки 

0,293475 

0,62907 

3,618187 

-0,181239 

Кустарники 

-1,202167 

1,113413 

2,90909 

-2,136832 

Вода 

-5,148586 

-4,03138 

4,551272 

-5,124147 

Мшистые 

болота 

-0,70012 

-0,823881 

3,101082 

0,596449 
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КИС, рассчитанный по значениям соответ- 
ствующих матриц ошибок составил 65 % для пред- 
ложенного подхода и 55 % для результата, получен- 
ного с помощью модуля СА_Магкоѵ. 

По результатам пространственного анализа те- 
матической карты на 2000 г., можно сделать вывод, 
что пойменные участки исследуемого фрагмента 
земной поверхности занимали общую площадь 
около 257,6 км 2 , а к 2002 г. по данным сформиро- 
ванной прогнозной карты - около 205,9 км 2 , то 
есть общая площадь пойменных участков умень- 
шилась приблизительно на 51,6 км 2 . В то же время 
мшистые болота на 2000 г. занимали общую пло- 
щадь того же фрагмента земной поверхности около 
1198,5 км 2 , а к 2002 г. - около 1258,1 км 2 , то есть 
можно говорить об увеличении площади данного 
типа земной поверхности на 59,6 км 2 . Водные 
объекты исследуемой земной поверхности, а также 
границы кустарников практически не претерпели 
никаких изменений, что хорошо видно из табл. 2 и 
3 и прогнозных карт (рис. 4). 

Заключение 

Предлагаемый подход к моделированию изме- 
нений земной поверхности, реализованный в ПО 
лаборатории ГИС ТПУ, использует в своей основе 
вероятностное формирование правил функциони- 
рования КА. Предварительные результаты иссле- 
дований, полученные с использованием модель- 
ных разновременных изображений, показали перс- 
пективность предложенного подхода. Результаты 
апробации при решении задачи анализа динамики 
земной поверхности фрагмента территории 
ХМАО, подтвердили результаты исследований, по- 
лученные на модельных изображениях, и показали 
основные тенденции изменений ландшафтных ти- 
пов изучаемой территории. 
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